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Assuming a large - scale homogeneous magnetic field, we follow the covariant and gauge - invariant 
approach used by Tsagas and Barrow to describe the evolution of density and magnetic field inho- 
mogeneities and curvature perturbations in a matter - radiation universe. We use a two parameter 
approximation scheme to linearize their exact non - linear general - relativistic equations for mag- 
neto - hydrodynamic evolution. Using a two -fluid approach we set up the governing equations as 
a fourth order autonomous dynamical system. Analysis of the equilibrium points for the radiation 
dominated era lead to solutions similar to the super - horizon modes found analytically by Tsagas 
and Maartens. We find that a study of the dynamical system in the dust - dominated era leads 
naturally to a magnetic critical length scale closely related to the Jeans Length. Depending on the 
size of wavelengths relative to this scale, these solutions show three distinct behaviours: large-scale 
stable growing modes, intermediate decaying modes, and small-scale damped oscillatory solutions. 



I. INTRODUCTION 

Magnetic fields play an important role in our Uni- 
verse. They appear on all scales, from the solar system, 
through interstellar and extra- galactic scales, to intra - 
cluster scales of several Mpc. Although magnetic field 
inhomogeneities have not yet been observed on scales 
as large as those exhibited by Cosmic Microwave Back- 
ground anisotropies, it is natural to expect that magnetic 
fields exits on such scales j|] , and that they could play a 
role in the formation of large-scale structure. 

Magnetic fields in intergalactic space are notoriously 
invisible. In the early 1950's the discovery of synchrotron 
radiation from the interstellar medium of the Milky Way 
led to the realization that it possessed a magnetic field. 
In 1957 Bolton and Wild |^ suggested that the Zeeman 
splitting of a radio transition should be observable in the 
interstellar gas. This provided a direct way of measuring 
the strength of this uniform magnetic field. That almost 
ten years passed between this suggestion and the first 
detections of this effect is an indication of the techni- 
cal difficulty such measurements posed at the time 
Since then, advances in observational techniques have 
produced firmer estimates of magnetic field strengths in 
interstellar and intergalactic space. 

We now know that a large scale, organized magnetic 
field fills the disk of the Milky Way. Studies of external 
galaxies indicate that all disc galaxies are permeated by 
large scale magnetic fields and that level fields 
are common in spiral galaxy discs and halos. Strengths 
of magnetic fields in intergalactic gas in 'normal' galaxy 



clusters have been measured using Faraday Rotation 
Measures (RM) combined with X-ray data. This gives 
typical field strengths for these fields of between 2 and 
G/iG, 1^, which is comparable to the field strengths in the 
denser interstellar medium in our galaxy. In 1993, Taylor 
and Perley Q discovered that these field strengths were 
exceeded by those in some less common, but more dense 
cooling flow clusters, where the field strength can be as 
high as 30fiG. These fields show ordered components on 
super - galactic scales |^ . 

Much of our local universe consists of voids, regions 
containing very little baryonic matter. As yet, the weaker 
fields within these voids, which may be relics of true cos- 
mic primordial fields, remain unmeasured. Primordial 
and protogalactic magnetic fields represent the large- 
scale fields which could play a role in structure formation. 
There are at present over thirty theories about the ori- 
gin of cosmic magnetic fields at galactic and intergalactic 
scales. Battaner & Lesch [Q look at astrophysical argu- 
ments to examine these models. These can be divided 
into four main categories, based on when the fields are 
generated, namely: a) during inflation, b) in a phase 
transition after inflation, c) during the radiation domi- 
nated era, and d) after recombination. It is the large 
scale fields which were produced during inflation which 
are most likely to have implications for structure forma- 
tion. 

Investigating the effect of magnetic fields on struc- 
ture formation is not a recent endeavour. In 1971, Ruz- 
maikin & Ruzmaikina Q gave a Newtonian analysis of 
the growth of density perturbations in a perfectly con- 
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ducting medium with magnetic fields. Wasserman ||] 
considered magnetic influence on galaxy formation, and 
angular momentum. Kim et al ^ extend this by includ- 
ing the back - reaction of the fluid on to the field, and as 
a result, derived a magnetic Jeans length. Battaner et al. 
Ip^ present relativistic analysis of the evolution of mag- 
netic fields and their influence on density inhomogeneities 
in a radiation dominated universe, while Jedamzik et al. 
and Subramanian & Barrow pl| have considered mag- 
netic dissipative effects at recombination. 

More recently Tsagas and Barrow p^-p^ have devel- 
oped a covariant and gauge - invariant approach to the 
analysis of magnetized density perturbations, in a uni- 
verse model containing a perfect fluid. Using the co- 
variant gauge - invariant formalism of Ellis and Bruni 
[ p^ , they derived a set of exact non - linear equations 
for general -relativistic magneto - hydrodynamic evolu- 
tion. On linearizing these equations around a flat Fried- 
man-Robertson-Walker (FRW) model, they were able 
to identify the general relativistic corrections to earlier 
Newtonian work, including the existence of a magneto - 
curvature coupling. 

Tsagas and Maartens extended these results, con- 
sidering also shape distortion effects due to the fields as 
well as density and rotational perturbations They found 
explicit solutions to the perturbation equations for the ra- 
diation and dust eras, as well as pure - magnetic density 
perturbations. They also identified and analysed other 
sources of the magnetic effects. 

When studying the evolution of density perturbations, 
the dynamics of the perturbed universe can be repre- 
sented by an autonomous dynamical system, described 
by a set of coupled differential equations for the gauge - 
invariant perturbation variables ||l7|-[l9[|. This allows the 
stability behaviour of the model to be investigated with 
relative ease, with no need to solve for the background 
variables. Furthermore, stationary points of the dynam- 
ical system correspond to exact or approximate analytic 
solutions of the linearised perturbation equations, thus 
providing a useful tool for obtaining solutions in cosmo- 
logically interesting situations. 



This paper is outlined as follows. In section III, we set 



up the evolution equations of the density and magnetic 
field inhomogeneities, and the curvature perturbations. 
We discuss in detail the approximations and linearisa- 
tion procedure in sections IV & ^ and identify an ex- 



tra term in the linearized prop agat ion equation for the 
spatial curvature /C. In section [Vl|, we set up a five di- 
mensional autonomous dynamical system equivalent to 
the linearized propagation equations derived in section 



III -these equations give a detailed description of the 
full phase - space of solutions for perturbation dynamics 
in a magnetized dust - radiation universe. The analysis 
and discussion of this system, with particular emphasis 
on the three dimensional invariant sets representing the 
initial radiation and final dust dominated states are given 
in section VII D. Finally, conclusions are given in section 



VIII. 



II. DEFINING VARIABLES 
A. Spacetime Splitting 

For all of the following equations, we will use standard 
units defined by SttG = c — \. In cosmology, the av- 
erage velocity of matter at each spacetime event defines 
a unique 4 -velocity vector: u°- — with UaU°' — — 1. 
The fundamental fluid- flow lines defined by this vector 
field are a congruence of worldlines, carrying the funda- 
mental observers. The projection tensor hab = gab+UaUb, 
(with h°'ch^i, = h"'b, and /i°a = 3) projects into the local 
rest spaces of comoving observers. 

The covariant derivative of any tensor T°'^cd may be 
split into a time derivative along the fluid flow. 



1 cd — U VeJ cd 



(1) 



and a totally projected covariant derivative orthogonal 
to the fluid flow [Bol 



(2) 



with total projection on all free indices. It is worth noting 
that if has non - zero vorticity, V is not a proper 3 - 
dimensional covariant derivative. Angle brackets denote 
orthogonal projection of vectors t^^"^ — h°'bV'', and the 
projected symmetric trace -free (PSTF) part of tensors 

Following lH], we can define a covariant spatial diver- 
gence and curl that generalize the Newtonian operators 
to curved spacetime: 

divy = v°K, (divr)a - v^Tab , 

CUrlK = £afccVV^ CUrlTafa = ecd(aV"T,)'^ , (3) 

where Sabc = Vabcdu'^ is the projected covariant permuta- 
tion tensor. 



B. Kinematic Quantities 

The projection tensor, together with the velocity 4- 
vector, split the covariant derivative of Ua into irreducible 
basic kinematic quantities |l2l ||20||: 



VbUa = -Ubila + "^bUa 



-Ubila + ^Qhab + (Jab + ^^ab 



(4) 



where iib = W aUbu"- , 6 = VqU°, aab = {aUb) and ujab = 
V[f,Ua]. These are, respectively, the 4 - acceleration, the 
expansion scalar, the shear tensor and the vorticity ten- 
sor. The expansion scalar is used to introduce a represen- 
tative length scale a along the observers' worldline which 
in turn defines the Hubble expansion rate, ^ = |0 = H. 
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The evolution of the expansion scalar is described by 
Raychaudhuri's equation: 

e-A + 2{a^ - uj^) + ie^ + i(p + 3p + B^) = Q, (5) 

where A = Vila, cr^ = ^(7"-''<Tab and = ^iLi"-''ujab (see 
|§ for details). 

Local curvature is described by the Ricci tensor Rab 
while non-local tidal forces and gravitational radiation 
are described by the 'electric', Eab and 'magnetic', Hab 
parts of the Weyl tensor. These are a pair of symmetric, 
traceless and completely spacelike second order tensors, 
defined by: 

Eab = E(^ab) = Cacbd.u'^u'^ , (6) 
Hab — H(^ab) = ^SacdC'^'^be'U''' ■ (J) 



C. The Electromagnetic Field 

The electromagnetic field is represented by the 
Maxwell tensor Fab = F[ab] ■ As seen by an observer with 
4 -velocity Ua, the field tensor splits into an electric Ea 
and a magnetic Ba part (see ||2^] for details): 

Ea = Fabu'' = -FabU'' , (8) 



divB = 2uj''Ea , (14) 

and 

a^'^'Bb - e^^^UbBr^ - ^QB" - e^^^ubEc 

- cnvlE" = B<''> . (15) 

The 4 - current J" enters into equations ( [T2| ) and ( |T^ ) 
through the charge density q = —J°'Ua and the projected 
current j° = J"^"^. The vorticity vector uja is related to 
the vorticity tensor via ujab = £abc^'^ ■ 

B. Gauge - invariant Variables 

Following ||l^ and ||l2|, the key covariant and gauge - 
invariant variables describing inhomogeneity in the fluid 
and the magnetic field are the comoving spatial gradi- 
ents of the energy density p, the expansion Q and the 
field density B^: 

Aa = ^^, ^aVaO, Ba^VaB\ (16) 

P 

where B^ — BaB"^. Other basic variables which appear 
in the exact equations are the spatial gradient of the pres- 
sure, Yq = VaP, the anisotropic pressure generated by the 
magnetic field, Hab — \B^hab — BaBb, and the comoving 
spatial gradient of the field vector: Bab = aVbBa- 



Ba = ^£acdF'^'^ , (9) 

where EaU°- — 0, and BaU°- = 0. 

III. EQUATIONS 

A. Maxwell's Equations 

Using the above decomposition. Maxwell's equations 
can be written as follows: 

VbF"^ = r , (10) 

and 

y[cFab] = , (11) 

where J° is the 4 - current which generates the electro- 
magnetic field and obeys the conservation law Va J° = 0. 

These equations can be split into spatial and temporal 
parts using the notation defined above: 

dWE + 2uj''Ba = q , (12) 

.^^Eb^e^^^..Eb~l^E^^e^^^UbB. 

+ curlB" = + f , (13) 



C. Medium and Propagation Formulae 

The matter description under consideration is a mix- 
ture of dust and radiation which share the same 4- 
velocity u", interacting only through gravity. In general, 
the evolution equations for density perturbations couple 
to an entropy evolution equation through the equation of 
state p = p{p, s), where s is the entropy density. How- 
ever, since entropy perturbations are only important on 
very small scales, we will assume adiabatic perturbations 
at all times, thus ignoring any entropy contributions [ p2| . 

The total fluid equations take the same form as those 
for a perfect fluid of infinite conductivity. The infinite 
conductivity approximation allows the electric field to be 
omitted from Maxwell's Equations while spatial currents 
are preserved . 

In this case. Maxwell's equations (^2|) - ( |l5| ) generate 
three constraints ]l6[ |, 

^^Ba = \q , (17) 

cm\B''^e''^''BbUc+ f , (18) 

divB = , (19) 
and one propagation equation 

^<a> ^ ^abQ^ _^ e^^^BbUJc - feS" . (20) 
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D. Conservation Equations 

Since we are only considering adiabatic perturbations, 
we only need the total fluid equations. The energy and 
momentum density conservation equations are respec- 
tively pi: 



p + pii + w)e = 



(21) 



and 



Pil+W+ '^)Ua +Ya- -B[ab]B'' + il^'Uba = , (22) 



where w — ^. Ignoring entropy perturbations, the pres- 



sure can be written as a function of the energy density 
only, i.e. p = p{p). It follows that the relationship be- 
tween Ya and Aa is given by 

aYa = pclAa , (23) 

where = ^ = ^ is the speed of sound in the fluid. 



E. Exact Propagation Equations for the 
Inhomogeneity Variables 

The propagation equations for the various spatial gra- 
dients are given in |I3]. They are, firstly, the equation 
for the comoving fractional spatial gradient of the energy - 
density: 



A<a> = wOAa - (a^ + Lu''a)Ab - (1 + w) Ga 

29 ^ . . 2aQB^ . aO 

-^[ab]0 



3p 



ha 



(24) 



secondly, the equation that governs the spatial gradient 
of the expansion: 

26 

e<a> = Qa - {<y\ + u:\)Qb - \pAa + BbaB^' 

- 3B[ab]B'' + aTZiia + ^aii'^Uba + aAa 

- 2aVa(cr2 -Cj2) , 

where Aa = ha''A.b = VaA, A = V'^Mq and 



(25) 



with 



+ A - 3{a'^ - Lu^) 



K = p+^B^-\e + a^-. 



(26) 



(27) 



representing the 3 - Ricci scalar of the observers instan- 
taneous rest space. 



Finally, the equation describing the evolution of the 
orthogonal spatial gradient of the magnetic field is given 
by: 

a-^ha^hb^ia^BepJ = -Bae{a% + C^^b) + {cTa" + i^a^B^b 

- iBaQb + 2aBea^[aUb] 

- 2aBeUj''(^aUb) + aB'^Vbic^ae + t^ae) 

- \Qa{2BaUb + UaBb) 

+ au''Be{Oab + i^ab + jQ^afc) 

- aha''RepbqBPu'> . (28) 



IV. DETAILED DISCUSSION OF 
APPROXIMATIONS 



In the papers by Tsagas, Maartens and Barrow - 
Jl^ , the cosmic magnetic field is treated as a coherent 
test field propagating on the background. The large- 
scale magnetic field Ba is assumed to be too weak to de- 
stroy the spatial isotropy of the background spacetime, 
which is taken to be a flat Friedman- Robertson- Walker 
(FRW) model^. This is an acceptable physical approxi- 
mation when the magnetic field density is small compared 
to the energy density of the fluid i.e. B^/p << 1. This 
issue is discussed in [g^ , where Zel'dovich calculates that 
B'^ / p < 8 X 10^^ if the model is to be applicable during 
the whole of the radiation dominated era. Furthermore, 
the isotropy of the microwave background can be used 
to place strong limits on the magnitude of the magnetic 
field H-l^. 

Other authors have used alternative approaches to the 
problem of maintaining isotropy in the background in the 
presence of magnetic fields. Battaner et al. [Q show that 
a mean magnetic field is incompatible with a Robertson - 
Walker metric and therefore assume that there is no mean 
magnetic field on cosmological scales, i.e. < >= 0. 
However, they include the presence of magnetic fields 
in smaller cells, with field directions random on larger 
scales. Thus, although there is no mean magnetic field, 
the model includes an average magnetic energy density 
< B^ 0. Kim et al. |^ also follow this approach, 
considering that, at recombination, field directions are 
randomly oriented on scales smaller than the Hubble ra- 
dius. 

In order to maintain the coupling between mag- 
netic irregularities and energy density perturbations in 
a straightforward way, a perfectly conducting medium is 
introduced. Looking at the covariant form of Ohm's law 
||12[, we have. 



^Spatial flatness is necessary for the 3 -Ricci scalar to be 
gauge - invariant. 
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(29) 



where a represents the conductivity of the medium. Pro- 
jecting into the rest space of a fundamental observer 
yields 



J" 



(30) 



To have a vanishing electric field, while maintaining non- 
zero spatial currents (J^"^ ^ 0), the conductivity of the 
medium must be infinite, i.e. a —>■ oo. 

An alternative approach would be to assume a pure 
magnetic field with no electric field and no spatial cur- 
rents. This would reduce ( |l3|) to 



e^'^UbBc + curlS" = , 
which is equivalent to 

uibBc] + lB[cb] = . 



(31) 



(32) 



Although this appears to couple the acceleration to the 
field gradient, when this relation is inserted into (p^, all 
the magnetic terms cancel out. This effectively decouples 
the magnetic field from energy density inhomogeneities. 



V. APPROXIMATION SCHEME 



In order to linearize the equations (|24|), ( P5| ) & (28), 
we introduce two smallness parameters. The first, ei is 
used to measure the extent to which the gauge - invariant 
variables deviate from zero (their value in a flat FRW 
universe). The other parameter £2 is a measure of the 
Alfven speed B^/p. 



G 



l^afcl 
|VaPl 

pH 



< ei 



\Hab\ 

|Va9 
i?2 



(33) 



where \Eab\ = {EabE'^^'Y'^ , |Vap| = {^aP^^pf^ etc. 

Tsagas extends this definition to a magnetized uni- 
verse by arguing that closeness to a spatially flat FRW 
spacetime is maintained when additional restrictions are 
imposed as follows: 



in„ 



\VbBa[ 

H\BA 



< ti and 



\K\ 



(34) 



In an exact flat FRW spacetime, all quantities of order 
ei vanish identically. Thus, (Tab , ^ab, Ua, {A, Aa}, {Eab, 
Hab}, {Aq, Sq}, {Bab, Ba}, {ICa = K} are ah consid- 
ered to be first order in ei. 

Note that although the magnetic field vector Ba is con- 
sidered to be a zeroth order quantity, its magnitude must 
remain small so that it does not disturb the isotropy of 
the background. This ensures that the the anisotropic 
pressure generated by the magnetic field, Habi is negli- 
gible in the background, and therefore Hab may also be 
regarded as first order in ei. 

Linearization of the above equatioons is implemented 
as follows: All terms higher than first order in ci are 
dropped, as well as terms higher than first order in the 
Alfven parameter 62- Terms like 6162 are however kept. 

At the end of the calculation, terms first order in 62 are 
dropped relative to zero order terms in the coefficients of 
quantities that are first order in ei. This is permissi- 
ble since the magnetic field is very weak (B^ / p << 1). 
However, this must be done last, otherwise terms may be 
dropped relative to others which could later vanish. 



A. Zero -order quantities 

The energy density p, pressure p and expansion Q do 
not vanish in the background spacetime. The magnetic 
field, Ba, is treated as a small test field propagating on 
the background. It follows that these variables can re- 
garded as zeroth order in our approximation scheme. 



B. First -order quantities [in ei] 

In order for the metric in some region of spacetime 
U to be written in a perturbed FRW form, the follow- 
ing inequalities must hold for the smallness parameter ei 



VI. LINEAR EQUATIONS 

A. Conservation Equations 

Linearization leaves equation (|2^) unchanged^ but 

slightly modifies the momentum density equation (|22|): 

p{l+W+^)Ua + Ya--B[ab]B'' ^0 . (35) 

6p a ' 

This, together with (p3|), gives a useful expression for the 
acceleration vector: 

^— ^ (-B^abiB^ - clAa) . (36) 

1 + W+ 1^ I a \P / 



In deriving this, VaW and VaC^ are treated as first order in 
ei , as in M 
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It follows that the divergence of the acceleration A is 
given by ^ 



C. Scalar Equations 



A = 



w - 



2B± 

3 p , 



(37) 



B. Linearized Propagation Equations 

The linearized propagation equations for the key 
gauge - invariant quantities are obtained by dropping 
terms of 0(2) or higher in both linearisation parameters. 
Note that terms of mixed order (eie2) and first order in 
62 are kept to start with (but may be dropped later). 

In these equations, the expansion scalar Q is replaced 
by the equivalent expression involving the Hubble pa- 
rameter i.e. 9 = 3H. Also, to first order, all projected 
time derivatives of quantities that are first order in ei are 
equal to their normal time derivatives (e.g. A<a> = A^). 

It follows that the energy density propagation equation 
(p3) becomes: 



Aa=3wHAa-il+w)ea 



6H 



2aHB^ 



(38) 



The equation for the expansion gradient (^) involves the 
spatial gradient of the divergence of the acceleration Aa- 
This can be determined using expression (|37|), the identi- 
ties for commutations between spatial gradients and time 
derivatives given in Appendix and the spatial gradient 
of the 3 -curvature K, given by 



aVaK = 2pAa + aBa - ^HQa . 



(39) 



Using these expressions, the evolution of the expansion 
gradient (to first order in ei) is 



-V^Aa - 



1 



- w - 



2B± 
3 p 



2p{l 



2 B- 

3 p 



V'Ba 



6cl{l + w) 



+ - + p{^ + -^ + lf) 



2B' 



3p(i + ^ + !fj ip[} + ^ + ,p 

4^2 



P^a 



B' 



aHV ujab 



aBa 



2 



3p(i + ^ + !f 



(40) 



Focusing on the growth, or decay, of density inhomo- 
geneities, the vector field A^ contains more information 
than is necessary. We can extract the required informa- 
tion by considering a local decomposition of the spatial 
gradient of Ao first introduced in 



3. 



aVbAa = Wa 



(41) 



where Wab = ^[ab] represents rotations of the density 
gradient Aq, Yiab = ^{ab) ^ ^hab describes the variations 
of Aa associated with pancake - or cigar - like structures, 
and finally A = A^a = oV^Aq is related to spherically 
symmetric gravitational clumping of matter. It is this 
scalar variable that is important when examining struc- 
ture formation. 

We also need to consider the following complementary 
scalar variables jl^ : 



Z = aVGa , B = 



52 



JC 



(42) 



which represent spatial divergences in the expansion gra- 
dient, the energy density gradient of the magnetic field, 
and perturbations in the spatial curvature respectively. 

As before, we linearize the scalar propagation equa- 
tions by dropping terms of order e\ or e^i but retain all 
terms that are first order in 62. 

Equations describing the propagation of these scalars 
are: 

3HB'^ 

A = 3wHA - (1 + w)Z + — B 



2p 



HB' 



-JC 



2ctB' 



2HZ - -A 
2 



^2 



3 P 



HA 



(43) 



l + w + %Sl 



3 p 



^2 



V'-B 



2^2 



3p 1 + w + 



2 

3 P 



B'B--- 



B' 



V2A 



pA 



-HZ, 



(44) 



This decomposition in analogous to the decomposition of 
the first covariant derivative of the 4 -velocity: (W). 



6 



-Z 



1 + w 



w - 



2B± 
3 p 



'AMI 

3 p 



IC 



(45) 



Finally, the completely linearized propagation equations 
for B and /C are 



B 



3(1 + w 



and 



(1 



w) 



Wj 



(50) 



AHci 



AHB^ 



3p( 



A + 



2HB^ 



If) P(l 



2 5^ 

3 p 



(46) 



2B'^H 



(l + w) 3p{l + w) p{l + w) 



B 



(51) 



Equations (|49|) and (|50|) are identical to those in ||16[ , but 
as a result of slightly different linearisation equation 
( pl| ) has an extra term compared to the corresponding 
result in llBl. 



D. Final Linearized System 

Combining equations (|4^-^), and finally dropping 
terms of order £2 with respect to zero order quantities, we 
obtain a second - order differential equation for the scalar 
energy density perturbations: 



A = -(2 
+ c^V^A 



icl - 6w)HA + i(l - 6c2 + 8w- 3w^)pA 
2B^H 



+ -c 



p{l+w) 
f i(2-3c2 



c^A-i(l-3c2 
3w)/3c2/C . 



2w)pclB 



(47) 



ft is not immediately obvious that the term containing 
is negligible with respect to the other terms involving 
A. In order to clarify this, we can write in terms of w 
by noting 

.2 d{cl) . 



dw 

(1 + iv) 



H 



(48) 



In what follows we will consider the above equations for 
values of ui e [0, i], but more specifically, we will look 
at the dust w 

From equation (B8|) we see that w = 



w 



J — and that 
Thus, the term drops away during 
both these eras and is negligible at all other times. 

Without that term, the second order propagation equa- 
tion for A becomes 

A = -(2 + 3c^ - 6w)ff A + i(l - 6cl + 8u; - 3u;^)pA 



1(2- 



\{\~?>^,+2w)pclB^ 
3c^ + 3u')pc^/C . 



■V'B 



(49) 



VII. DYNAMICAL SYSTEMS ANALYSIS 

Instead of attempting to find exact solutions to equa- 
tions (49-^) for a dust - radiation background, we in- 
stead follow Bruni and Piotrowska and perform a 
qualitative analysis of the above system of differential 
equations. 



A. The System 

It is useful to change the independent variable from 
proper time to a function of the scale factor a: with the 
choice T = Ina (which yields dr/dt — H), and using the 
standard harmonic decomposition described in Appendix 
^, equations (p9|), ( [50| ) and ( ^ become, 



A" 



-i(l + 6c2 - l5w)A\^. + 1(1 - 6cl + 8w- 3w^)A(^ 



2 .2 p 



c^^4^A(,^-l{l-3cl + 2w)clB^,-, 



+ (2 - 3cl + 3w)c^/C(fc) 



-a;,, + :/ / A(fc) , 



^''^ 3(1 + w) f'^) (l + w) 



(52) 



(53) 



and 



4c2 

(rTt)^('=)"3(Tf^^(^-) + (TT^ (^^^ 

where a prime denotes differentiation with respect to r. 
From now on we will drop the subscript (fc). 



4c 



2c2 



''At this point Tsagas and Maartens [jl^l expand the LHS of 
equation ( [si] ) using the definition of IC (see equation (^)): 
K. = 2HIC + a?K and then drop the second term on the RHS 
of (|^) relative to 2HK.. 
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B. Coefficients in Terms of w 

In order to close the system, propagation equations are 
needed for all background variables. It turns out however 
that the coefficients of the above differential equations 
can all be written explicitly as functions of the equation 
of state parameter w, thus only an evolution equation for 
w needs to be found. To achieve this, we normalise the 
scale factor a at dust - radiation equi - density by intro- 
ducing the variable S = a/aE (see Ehlers & Rindler ||29| ] 
and Padmanabhan [Q). 

Since the two fluids are coupled only through grav- 
ity, the energy conservation equation is obeyed sepa- 
rately for each component. This gives pd = ^S~^ and 
Pr = ^S^^ for dust and radiation respectively, where 
is the total energy density at equi - density. Total en- 
ergy conservation yields p = ^{S~^ + S~^). The only 
contribution to the total pressure comes from radiation 
component: p — Pr — ^ — ^S~^. It is now straight- 
forward to find w in terms of S [n9| : 



1 



w 



3(5-1-1) 
which can be inverted to give 
1 - 3w 



S = 



(55) 



(56) 



In this way the expansion of the universe model is pa- 
rameterised by w, with w G [0, and varies from a pure 
radiation- dominated (t — > 0) to a pure dust - dominated 
{t — > oo) phase. 

These results ( |55| ) and ( [56| ) allow us to write the coef- 
ficients that occur in the perturbation equations as func- 
tions of w only (see Appendix |^). 

The energy density of the magnetic field has a ra- 
diation-like propagation equation = Be^S~^, hence 
the Alfven speed can be written as 



2 c 2 

c„ — bwc 



aE 



(57) 



The propagation equation for w can now be written 
down. It is 



w — 3w{w — -i) . 
In terms of w, (^) becomes 

A" = q^aA' + /3aA -t- 7aB + tjaJC 

where 



(58) 
(59) 



OA = -2- 



(1 



W) 



(1 



'A 



7A 



1 



w^{5 ~ 3w) 
(l + w) 
1 — w 



2w^ 



rjA = Gwc 



aE 



(2 



(l + w) 
- w + 'iw^) 



w) , 

8(l-3w)2 , 
9(1 + w) ' 



^E : 



(1 + w) 



(60) 



Similarly, ( |53| ) and (^ become 

B' = aB^.' + Pb^ , 

and 

K.' = a^A + PkJC -f -ixB , 

where 
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(61) 
(62) 



OLE 



Pb = 



UK 



3(1 + w) ' 

4w 
3(1 + w)2 

IQw 
3(1 -f w)2 
8w 



(1 - iw) , 



K 



IK 



^2 



(1 + w) " 

12w 2 

:^aE ; 



(63) 



{l + w) 

with the evolution of the background determined by (|5^ 

C. Reducing the Order of the System 

In the simplest case of a single fluid, equations (|5S|), 
( |6l| ) and ( |6^ ) together with the condition w = constant, 
form a fourth - order autonomous system of differential 
equations. In the more general two - component case, the 
system is five dimensional due to the inclusion of the 
propagation equation for w (58). However, the order of 
the system can be reduced by introducing a new vari- 
able U ll^. This reduction makes it possible to easily 
determine the qualitative behaviour of A, rather than its 
exact evolution law. We start by defining AT = A' and 



U = A/A = A'/A , 
n = VA2 + A2 . 



(64) 



7?. represents an 'amplitude' of the perturbation and it 
is not directly relevant to our analysis. In order to keep 
the system dimensionally consistent, we define two new 
variables, 



V ^ — , and W ^ — . 
A A 



(65) 



In terms of U,V,W, the equations (59), (^ and (|6 
become 

U' ^ -U^ + uaU + Pa + IaV + vaW , (66) 



V = asU + I3b-VU , (67) 



W ^ UK + PkW + -ikV ~WU , (68) 



with the coefficients determined by (|6C|), fc (|63|) , together 
with the propagation equation for w 
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1 (J ii 1 III/ 


VnhtP<i (in TI V W\ 


Eigenvalues (\i,2,3,a) 


Rl 


(io,-^,-^) 


1 

3 




2 


-l + 2c^B 


R2 


aE aE 

(i,-4cL, 1,-1(^ + 3)) 


1 
3 


4cL 


2 + 8c^E 


-l + 2cL 


R3 
R4 


\,-1 + 2cIeA,-{'^+7cIe) 
(i,2,l,i(l + 2cL)) 


1 

! 

3 


1 - -iclE 

-2 


1 - 6C^E 
-2 - 4cL 


3 — 2C„£; 

-3 + 2cIe 



TABLE I. Table of equilibrium points for 111 = ^ 



D. Analysis 

The system is non- linear, however we can analyse it 
locally by linearising about any stationary points, with- 
out losing the details of the qualitative dynamics |^,^ . 

From equation (p8) we can see that stationary points 



exist for w = ^ and w — Q. 

Since equation ( p8| ) decouples from the rest of the 
system, the 3 - dimensional subsystems corresponding to 
these values of w correspond to invariant sets describing 
the early radiation - dominated and late dust - dominated 
periods of dynamical evolution. We now consider these 
cases separately. 



E. Radiation Era: 



On substituting w = ^ into (|6y) and 
for the radiation- dominated era becomes: 



), the system 



U 



W = 1- 2cIj,W + Zclr,V -WU 



(69) 



The equilibrium points for ( |69| ) can easily be determined 
and are given in table | above. 

The 3 - dimensional space described by {U, V, W}, and 
w = ^ (see figure I) is an invariant set, so we can classify 
the equilibrium points according to the eigenvalues A2,3.4 
(since Ai = ^ is not relevant). 

Point Rl is a saddle (not included in figure I because 
it is far away from the other points). Orbits close to it 
may initially evolve towards it, but end up evolving away 
again. At Rl U — 0, so density perturbations neither 
grow or decay. Using equations (|6J) and ( |65|) it follows 
that A, B and K. are constant: 

A = const = Co , 



(70) 



Point R2 is also a saddle, but this time since U < 0, 
density inhomogeneities are decreasing. Solutions at R2 
again follow from equations (|6l) and (Bsf): 



A : 

B 

K, ; 



-4d 



2c- 



1 n ( a 



(71) 



Point R3 is a node source, representing unstable equi 
librium, so all orbits close to this point evolve away from 
it. Again since U — —1 
represent a decreasing density inhomogeneity 



2c^£ < 0, solutions at R3 



A 
B 

K. 



A = C2(^)(-^+'^-) , 

-(l + 7c2^)C2(^)(-^+'^-) 



(72) 



Point R4 is a stable node, or a sink. Orbits close to this 
point evolve towards it as (^) increases. Since U > so- 
lutions at R4 represent growing density inhomogeneities: 

A = C3(^)\ 

/C = i(H-2c^^)C3(^)^ 

(73) 

Solutions at the points given above correspond to ap- 
proximate solutions (to leading order in c^^) of the per- 
turbation equations during the radiation - dominated era. 
Therefore by linearity, the general solutions for the per- 
turbation variables A, B and /C are given by a linear 
combination of (^- ^3|) : 



A^Co + Ci(^)-^^-+C2(^)^ 
B = -^^Co + C^{^)-'^'^ 

+ C2(^)(-^+^^-) + C3(^)\ 
^--i^o-5^(^)-^^- 
-(l + 7c^^)C2(^)(-^+'^-) 

This can be written more concisely as: 



(74) 



9 



Point 


Values {w, U, V, W) 


Eigenvalues (Xi 234) 


Dl 

D2 
D3 


(0,0, ^[-^-|], const) 

(0,[-i+?],|,0) 
(0,[-i-C],|,0) 





TABLE II. Table of equilibrium points for w = 0, with ^ = iW6(i| - p^). 



w 




The equilibrium points of this system are shown in table 
||. In this case, the U and V propagation equations de- 
couple from the equation for W' , which means that the 
system is effectively only 2 - dimensional. A critical scale 
Xec = '2'^0,/kEc appears in the analysis of these equation 
through its corresponding wave number kEc- 



-130 



Figure 1: 

a 

a 



C(a) 



where a solves the cubic equation: 



f a Y 



(75) 



{2ci 



I)a2_2«-8c2s(l-|-c2^) = 0, (76) 



and corresponds to the super - horizon solutions given in 
, with slightly modified exponents - due to the extra 
term in the equation for the spatial curvature IC (p4). 



F. Dust Era 

When solving for the equilibrium points in the dust- 
dominated era, we substitute w — Q into the equations 
defining the dynamical system (p6^-pS). This yields: 



U' = -U^ - \U 
V = ^U-VU , 
W' = -WU . 



3 8 ;„ ^ 

2 g'^-E 



2<EkE'V , 



(77) 



27 



1 



16 1 + 3cIe 



27 
16 



(1 - 5cIe) 



(78) 



The solutions at each point will have different behaviour 
depending on the value of ks- In what follows, we first 
look at the general properties of the equilibrium points 
and their corresponding solutions and then give specific 
information about their behaviour for the three regimes 
of fcfi values that emerge from the analysis. 

Point Dl: These solutions are confined to the {U, V} 
plane for constant W and approach a constant value in- 
dependent of the value of fc^;. 



A 
B 



const = Co , 



2c^ 



l-2fe£ 



IC ~ const — Ck 



(79) 



Point D2: This point has different stability behaviour 
for different values of fcs . The curvature variable W van- 
ishes independently of fc^. Since the solution at the point 
depends on ^ = ^(fc_B), the nature of the solution mode 
(e.g. growing/decaying/oscillation) will depend on the 
value of kE- 



A = Cl(^ 

/c = o. 



(80) 



Point D3: This point has different stability behaviour 
in two different regions of ks- In both regions, the cur- 
vature variable vanishes. The solutions at D3 are: 



A = C^2(^) 



6^ 4^2(J^)(-4-€) , 

/c = o. 



(81) 
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As in the radiation dominated era (discussed in section 
VII E ), solutions (to leading order in c^^) for the per- 
turbation variables A, B and IC are given by a linear 
combination of the solutions at Dl, D2 and D3: 



B 
IC 



]Co + |Ci(-^) + |C2(-f^)-2 



(83) 



A = Co+Ci(^)(-4+«)+C2(^) 



In the dust era, ( — ) = (t^)^^'^, so we can write the above 
solutions ( p3| ) in terms of {j^) as in For example, 

,the solution for density perturbations is: 



/C = const = Ck 



|]^^0 + |Ci(^)(-4+«) 



^2(^)(-^- 



(82) 



2/3 



C2(f) 



We now look more closely at how the nature of the sta- 
tionary points Dl, D2 and D3 and their corresponding 
solutions depends on the wavcnumbcr ks- 

Region 1 - ks < ksc- See figure II 

Point Dl: In this region, the eigenvalues are, respec- 
tively, A2 = 0, A3 > 0, and A4 < 0. It follows that Dl 
corresponds to a line of saddles in the solution space for 
{U,V,W}. 

Point D2: Here, all the eigenvalues are either neg- 
ative or zero, so this point is a stable sink. Since 
^ > i ^ ^ > 0, it follows that the solution at D2 
corresponds to a growing density inhomogeneity. 

Point D3: All eigenvalues are positive, so this point 
is an unstable node. Since U = —j — < 0, it follows 
that the solution at D3 corresponds to a decaying density 
inhomogeneity. 



W 




It consists of a magnetic field induced constant mode, 
plus the two usual non - magnetized adiabatic modes. 
Comparing to [l6| equation (65), we see that the dynam- 
ical systems analysis does not recover the non - adiabatic 
decaying mode. This is because our analysis is set up to 
look for asymptotic solutions, so as i — > 00, the magnetic 
energy density decays faster than that of the fluid, and 
0. This removes the magneto - curvature coupling 
from the propagation equation for IC, and therefore our 
analysis does not obtain the extra decaying mode in- 
duced by the magnetic field. 

Region 2 -ksc < ks < 2^^^'^' ^"^^ figure III 

Point Dl: In this region, the eigenvalues are, respec- 
tively, A2 = 0, A3 < and A4 < 0. It follows that Dl 
represents a line of sinks (only one is shown in figure III) 
in the solution space for {U, V, W}. 

Point D2: In this region, the eigenvalues are, respec- 
tively A2 < 0, A3 > and A4 > 0, so D2 is a saddle point 
(or a line of sources if A2 = 0) in the solution space for 
{U,V,W}. It follows, since < 0, the solution at D2 
represents a decaying inhomogeneity. 

Point D3: Here we obtain the same behaviour as in 
the case when < kEC, i-e. D3 is a saddle, correspond- 
ing to a decaying density inhomogeneity. 



U 



Figure II: 

For very long wavelength solutions i.e. ks << ksc we 
have ^ ~ I . In this case the explicit solutions of the per- 
turbation equations (which are a linear combination of 
the solutions at the points Dl, D2 and D3) for the dust 
era are given by: 

A = Co + Ci(^) + C2(^)-i 



0.5^ -0.35_ -0.2 -0.05 O.I 
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Figure III: 

As for region 1 we can now write down explicit solu- 
tions of the perturbation equations. In terms of proper 
time -r- we obtain: 



Co + Ci(^)(-6+«')+C2(7^)(-6-«') 



" 2^ 



r 3 



-|]C0 + |Ci(^)(-6+«*) + |C2(^)(-6 



/C = const = Ch 



where 



6 




'^EC 



(85) 



The spiral behaviour round the stable point is evident 
from the damped oscillatory solution (see figures IV and 
V). 

In this region, the wavelength of solutions falls below 
a critical wavelength which is related to the Jeans length 
with a correction due to the magnetic field (see section 
VII G , below). Thus, these solutions do not result in 
^it) prowing density inhomogeneities, but in general, oscil- 
late as sound waves. The constant density solution due 
to the magnetic field (corresponding to Dl) represents 
(84^he only stable solution in this region, and can be seen 
as a spiral point in figure V. Note that this is the only 
equilibrium point in this figure, as the other 2 solutions 
represent pure sound waves. 



and < < i 

In this region, none of the solutions correspond to 
growing modes for density inhomogeneities. The line of 
stable sinks represent modes where the density pertur- 
bation approaches constant value. 

Region 3 -kE > i^/^^Ec- See figures IV and V 

Point Dl: For these values of fc^, the eigenvalues A3 
and A4 are imaginary. Since their real parts are both neg- 
ative, this solution acts as an attractor, and gives rise to 
a line of stable spiral points. This repetition of the equi- 
librium points is, as mentioned before, a consequence of 
one of the eigenvalues vanishing. Here, the density and 
magnetic parts of the solution decouple from the curva- 
ture {W equation), and thus the spiral point exists in 
any W — const plane (see figure V). 

Point D2: The solutions for U are imaginary in this 
region, so density inhomogeneities oscillates as a sound 
wave, neither growing nor decaying. 

Point D3: Here there are no real solutions. Density 
perturbations again oscillate as sound waves. 

This time the solutions are given by 



A = Co + Ci(^)(-4+'^) + C2(^)(-3- 



i/3) 



B 



1 



IC^Ck. (86) 



with /5 = 1 



we obtain: 



6(1 



kl 



In terms of the proper time 



+ (^^)-6 (^3C0s[|/31n(^)] + C4sin[|/31n(^)]) 



B 



1 



+ |(u)"MC3 cos[|/3 ln(4)] + C4 sin[|/3 ln(^)]) , 
lC = Ck. (87) 




2 

W 4 
6 




Figure IV: 



Dl 



U 



-0.1 -0.06 -0.02 0.02 0.06 0.1 



10 
7.5 
W 5 
2.5 
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Figure V: 



G. Magnetised Jeans length 

It is clears from the above discussions that the critical 
scale ksc which appears in our analysis is not quite the 
scale which determines the onset of oscillatory behaviour. 
This is determined by the wavenumber kosc = '^^^ec, 



whose corresponding wavelength is: 



3V3^ 

16727r 1 /-I 
15 He^ 



1/21 



3^2 



This scale is closely related to the Jeans Length as it de- 
termines the wavelength which divides oscillatory from 
growing or decaying solutions. It takes the form of a gen- 
eral critical wavelength, modified by a linear factor due 
to the magnetic field. It is clear that since the magnetic 
field contributes to the overall pressure, it enlarges the 
size of regions able to resist gravitational collapse, and so 
increases the critical scale relative to the non - magnetised 
case. 

In the limiting case c^^ 0, we obtain: 



An 



16V21T 1 
15 Hp, 



(89) 



which is the corresponding scale for a flat dust - radiation 
model (with no magnetic field) found by Bruni and Pi- 
otrowska in |jl^. 

The scale Xosc found above (^8|) differs considerably 
from that found in [T^ ] in their analysis of a magnetized 
dust solutions (see equation 68). In a pure dust model, 
the magnetic field acts as an effective pressure, allow- 
ing for the extraction of a 'magnetic Jeans length'. This 
length is directly proportional to the Alfven speed, and 
vanishes if there is no magnetic field. Furthermore it is 
considerably smaller (of order c^^) compared to the crit- 
ical length found by us using a two -fluid model, since 
in this case, the magnetic field is not the only source 
of pressure; there is a contribution due to the radiation 
component. 



VIII. SUMMARY 

Assuming a large-scale homogeneous magnetic field, 
we use the exact general relativistic propagation equa- 
tions given in jlj] to describe the evolution of density 
and magnetic field inhomogeneities and curvature per- 
turbations. We use the usual approximations of small 
magnetic field energy density {B"^ / p << 1), and infinite 
conductivity (cr — ^ oo E' — > 0), to simplify these equa- 
tions. We carefully linearise the exact equations using 
a two - parameter approximation scheme and find an ex- 
tra linear term (compared to 13 1) in the propagation 



equation for the spatial curvature IC. This extra term 
however makes no qualitative difference to the dynamics. 

Rather than attempting to solve these equations ana- 
lytically, we follow Bruni and Piotrkowksa jl9| by set- 
ting up the linear perturbation equations as a four- 
dimcnsional autonomous dynamical system. This ap- 
proach provides an elegant description of the dynam- 
ics of gravitational and magnetic inhomogeneities for a 
two component medium comprising of pressure - free mat- 
ter (dust) and radiation, interacting only through grav- 
ity. The equilibrium points for the three - dimensional 
invariant sets representing the radiation (w = ^) and 
dust {w — 0) eras correspond to the approximate solu- 
tions of the perturbation equations derived in Tsagas and 
Maartens Jl^l. For example in the radiation era we find 
the following solution for density perturbations on scales 
much larger than the Hubble radius: 



(-l+2cfB) 



+ C,{f-f . (90) 



The adiabatic decaying mode decays even less rapidly due 
to the magnetic field, and the decaying non - adiabatic 
mode decays slightly faster than the corresponding so- 
lution in |n6[ . To zeroth order in c^, these solutions 
correspond exactly to those given in pa (see equations 
(50) and (54) in their paper). 

The extra term found in the linearisation of the cur- 
vature propagation equation is proportional to c^^, and 
thus tends asymptotically to zero. Since the dynami- 
cal systems analysis examines the asymptotic behaviour 
of these models, this extra term makes no difference in 
the analysis of the dust era. Similarly, the superhorion 
dust solutions do not include the non - adiabatic decay- 
ing mode induced by the magnetic field, since the Alfven 
speed (c^ — j p) also tends to zero at late times. 

An important feature that arises out of the dynamical 
systems analysis, is that we obtain three distinct evolu- 
tion regimes for the perturbation modes. These regimes 
are defined through a critical scale A ec which is of the 
order of the Hubble scale at matter - radiation equiden- 
sity and includes a modification linear in c^^ due to the 
presence of the background magnetic field: 



A 



■EC' - He 3^3 



The three regimes are: 

• Large wavelengths -\e > ^ec These are large- 
scale density and magnetic inhomogeneities that 
grow unbounded giving gravitational instability. 
Solutions for these wavelengths consist of a grow- 
ing mode, a decaying mode, and a constant mode, 
represented by a saddle point in the phase space of 
solutions. 

• Intermediate wavelengths - Xosc < Xe < X ec 
These perturbations are over - damped and there- 
fore decay asymptotically to a constant value Cq 
without oscillating. 
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• Small wavelengths -A^; < Xosc Here the perturba- 
tions oscillate like sound waves, while their ampli- 
tude decays. 

Furthermore, the analysis of the dust models yields a 
magnetic field corrected scale Xosc closely related to the 
Jeans length, which is a multiple of the general critical 
scale Xec for flat radiation - dust models: 



(1 + He) , 



15 He 



(91) 



where the quantities c^^. He, Xe are evaluated at equi- 
density. 

Finally we note that Xosc is more general than the 
'magnetised Jeans Length' found in [|l6[ for a pure dust 
model, since it takes into account the pressure effects re- 
sulting from a proper two - fluid description. 



APPENDIX B: HARMONIC DECOMPOSITION 



When setting up and analysing the dynamical system, 
all the equations can be reduced to ordinary differential 
equations if we restrict our attention to the harmonic 
components of the perturbation variables. This is a way 
of effectively separating the time from the space vari- 
ables, and involves writing the perturbation scalars in 
terms of harmonic scalars as follows [H : 



A = ^ A(„)Q(„) , 

n 

S = ^ B(n)Q(n) , 
n 

/C = ^ K,{n)Q(n) ■ 



(Bl) 
(B2) 
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The scalar harmonics, Q, are defined by 

Q(n) =0 and 



(B3) 
(B4) 



APPENDIX A: COVARIANT IDENTITIES 

This identities are written using the notation of ||2l| ] 
and ||l^, except using V instead of D for the orthogo- 
nally projected covariant derivative. These identities are 
used in deriving the propagation equations (assuming a 
fiat background and vanishing cosmological constant): 



curlV„/ 

V2(Va/) 
(aVaJfo..)- 

div curl V 
(div curl S)a 
curl curl Va 
curl curl Sab 



-2fuja , (Al) 

ay J + afAa , (A2) 

V,(V2/) + 2/curla;a , (A3) 

aVaife... , (A4) 

= V[,V^,]^=^ (A5) 

, (A6) 

icurl(div5), , (A7) 

Va(divF)-VVa (A8) 

|V(,(div5)b) - V^Sah , (A9) 



where the vectors and tensors vanish in the background 
and Sab = S(^ab) ■ The magnetic field vector does not van- 
ish in the background, and so its projected derivatives do 
not commute to li near order. For the magnetic field, the 
vector identity in (A5) must therefore be changed to: 



V[aVb]i3c — \'R'dcbaB'^ — Sabd^'^Bc , 



(AlO) 



where TZabcd is the 3 - Curvature tensor formed from Rabcd 
and the kinematic quantities |l^,p^. 



where n = fc > (since we are dealing with a fiat model). 
Spatial fiatness of the background also means that k is 
simply related to the wavelength A of the perturbation 



since A = 



APPENDIX C: COEFFICIENTS IN TERMS OF W 



cj. is only formally the speed of sound since the fiuids 
interact only through gravity. Using the energy - density 
conservation equations for matter and radiation we can 
express in terms of the parameter w: 



4w 



p 3(H-w) 



(CI) 



The total energy density of the matter - radiation mixture 
is given by: 



27 u>3 



P = -K-PE 



2'^ [I- ' 



(C2) 



where pE is the energy density at matter - radiation 
equality. This leads to a similar equation for the Hubble 
parameter H^: 



= -:rH 



27 „o 



2 ^ (1 - 3m;)4 



(C3) 
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where He is the value of the Hubble parameter at mat- 
ter - radiation equality. 

When harmonically decomposing the perturbation 

equations in section 



VII 



the Laplacian terms give rise 
to coefficients of the form /a^H^. These can be ex- 
pressed in terms of the variable w as follows: 



fc2 2(1 - 3m;)2 



2^2 



E 

2(1 



3w 



3w 



(C4) 



where k^- 
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